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The effects of image charges (i.e., induced surface charges of polarization) in spherical geome- 
try and their implication for charged colloidal systems are investigated. We study analytically and 
exactly a single microion interacting with a dielectric sphere and discuss the similarities and discrep- 
ancies with the case of a planar interface. By means of extensive Monte Carlo (MC) simulations, 
we study within the framework of the primitive model the effects of image charges on the structure 
of the electrical double layer. Salt-free environment as well as salty solutions are considered. A 
remarkable finding of this study is that the position of the maximum in the counterion density (ap- 
pearing at moderately surface charge density) remains quasi-identical, regardless of the counterion 
valence and the salt content, to that obtained within the sin(;/e-counterion system. 
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I. INTRODUCTION 



II. THEORY 



In charged colloidal systems electrostatic effects, and 
especially the structure of the electrical double layer, 
often play a crucial role in determining their physico- 
chemical properties. It is well known that charged col- 
loids (i.e., macroions) have typically a low dielectric con- 
stant (£r ~ 2 — 5) which is much smaller than that of the 
surrounding solvent (e.g., for water ~ 80). In most 
of the theoretical works, this dielectric discontinuity is 
ignored. 

Nevertheless, a few studies have addressed the effects 
of image charges (i.e., image forces stemming from the di- 
electric discontinuity) on the counterion distribution for 
planar geometry which is closely relatejito our problem. 
An electrolyte close to a charged wallldu or confined be- 
tween two charged platesa had been the subject of MC 
simulations. Similar|-swS|tems have also been invesljiffated 
by integral-equationQBEl and mean field theories .QotI 

As far as the spherical geometry is concerned, much 
less literature is available. Counterion distributions with 
image forces in salt-free environment had been investi- 
gated by MC simulations .1121 There an oversimplified ap- 
proximation for the treatment of the image forces was 
used. The main conclusions however remain correct on a 
very qualitative level of description. 

The aim of this paper is to provide a detailed analysis 
of the image forces in spherical geometry and their ef- 
fects on the structure of the electrical double layer. The 
remainder of this article is set out as follows. Section II 
corresponds to the analytical part of the paper. We first 
briefiy present the general theoretical background of the 
concept of image charges in spherical geometry. Then we 
apply it to colloidal systems to compute (exactly) some 
relevant observables and discuss our results. Section III 
is devoted to the computational details of our MC sim- 
ulations. In Sec. IV we present our simulation results 
for salt-free environment as well as salty solutions where 
image forces are explicitly taken into account with no ap- 
proximation. Finally, Sec. ^ contains brief concluding 
remarks. 



In this part we mainly study the interaction of a single 
excess charge with a dielectric sphere. We briefly present 
the formalism of the dielectric model for spherical inter- 
faces and discuss some important electrostatic proper- 
ties. Such a system captures the underlying physics of 
image forces in spherical geometry. Moreover, a system- 
atic quantitative comparison with the planar geometry is 
undertaken. 



A. Poisson equation with azimuthal symmetry 

The model system is sketched in Fig. |^. Consider an 
uncharged dielectric sphere of radius a and dielectric con- 
stant (relative permittivity) £2 embedded in an infinite 
dielectric medium (region 1) characterized by £1. A sin- 
gle excess charge of magnitude q is located outside the 
dielectric sphere at a distance h = |b| from its center. 

The central problem is to determine the electrostatic 
potential ^(r) at any point in the space. This is achieved 
by solving the Poisson equation which reads 



A$(r) = - 



P(r) 



(1) 



where p(r) is the volume charge density and e = eoSi with 
£0 being the vacuum permittivity and i — 1,2. Since here 
p{r) ~ qS{r — b) and taking into account the azimuthal 
symmetry, Eq. (0) reduces (for r ^ b) to the Laplace 
equation 



A$(r, ( 



J.2 Qj, 
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where 9 is the angle between r and b (see Fig. |i|) and 
r = |r|. The general solution of thp-|^pl^ce equation 
with azimuthal symmetry is given by 
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FIG. 1: Model for a dielectric sphere (colloid) of dielectric 
constant £2 embedded in an infinite medium characterized 
by a different dielectric constant £1. An excess charge (q) is 
located near the boundary outside the spherical particle. This 
is a two-dimensional representation of the three-dimensional 
system. 
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Piicos0), (3) 



where Pi (cos 9) is the associated Legendre polynomial of 
order /. 

Inside the dielectric sphere (region 2) the electrostatic 
potential $2(1") must be finite at r = so that A^; = in 
Eq. (p|), and hence 



B. Boundary conditions 

The electrostatic potentials given by Eqs. (|^) and (|^) 
will be univocally determined by applying the proper 
boundary conditions that will fix Ai and Ci . The bound- 
ary conditions are derived from the full set of Maxwell 
equations. The results are that the normal components 
of the displacement D and the tangential components of 
E on either side of the spherical interface at r = a satisfy 



(D1-D2) -1112=0 

(El - E2) X ni2 = 



(7) 



where ni2 = r/r is a unit normal vector to the surface 
directed from region 2 to region 1 (see Fig. |l]). Within the 
framework of the linear response theory we have D = eE. 
Combining Eqs. (^ and (||) with Eq. (^ and noting that 
E = -V$, it follows that 
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This set of two equations [Eq. (||)] can be readily solved 
to yield the Legendre coefficients Ai and C/: 



A, 



Ci = 
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and hence 
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Concerning the electrostatic potential outside the di- 
electric sphere (region 1) we know that without dielectric 
discontinuity (at r = a) the potential would simply be 



given by 

tity 



47reoei |r — b| " 



Making use of the following iden- 
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where r< (r>) is the smaller (larger) of r aad 6, the elec- 
trostatic potential $i(r) in region 1 readsM 



$i(r 
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recalling that $i(r) must be finite at r ^ 00 so that 
Mz = in Eq. (PI). 
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The physical interpretation of Eq. (^^ is straightfor- 
ward. The first term represents the usual electrostatic 
potential (without image forces) generated by q and the 
second term can be referred to as the electrostatic poten- 
tial due to "image charges" stemming from the dielectric 
discontinuity. As expected, the strength of the image 
force is strongly governed by the jump Ae in the dielec- 
tric constant defined as 



Ae 



£1 - £2- 



(11) 



In particular, one can anticipate and state that the inter- 
action between the microion q and the dielectric particle 
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(i.e., the self-image interaction) is repulsive for Ae > 
(i.e., El > £2) and attractive for Ae < (i.e., £1 < £2) as 
it is also the case in planar geometry. 

One can show th£it,Eq. (^o|) can also be written as 
follows (see e.g., Ref.ll3 and references therein) 



$i(r. 
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where u = ha^/b^ (see Fig. pjlj In this formalism the 
geometrical structure of the image charges is transpar- 
ent and it is specified by the second main term (between 
brackets) of Eq. (p^). More precisely, one has to deal 
with an infinite manifold of image charges distributed 
along the oriented segment u that electrically compen- 
sates the image point-charge qim located at u and whose 
magnitude is given by 



The critical angle 6* where cr^^^'''' changes sign is given 
by the geometrical condition 



Ei(r = a,r) ±ni2 
E2(r = a,r) _Lni2 



(17) 



which is the orthogonality condition at the interface be- 
tween the (inner and outer) electric field and ni2. In 
terms of Legendre polynomials, Eq. ( pT|) can be equiva- 
lently written as 
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where Eq. ( jig ) was used. Two limiting cases can be 
easily described: (i) for b/a 1 we have 9* — > 7r/2 
[recalhng that Pi (cos 0) = cos^] and (ii) for 6/a ^ 1 we 
have 9* 0. In general, 9* increases with b and it is a 
complicated function of b/a, £1 and £2. 

For a planar interface, the surface density of polariza- 
tion charge (j'"^ot^\d) is given by 3 



£1 - £2 a 
£1 + £2 6 ' 



(13) 



C. Polarization charge 

It is important to know the surface distribution of the 
induced charge on the spherical interface. In the bulk 
(i. e., in region 1 or 2) we have a zero volume density of 
polarization charge [ppoi) since ppoi = £oV-E = — V-P = 
(except at r = b). At the interface (r = a) the surface 

density of polarization charge cr^^^''^ is given by 
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where d — \/ x^ +2/^ is the radial distance (in cylindri- 
cal coordinates system) belonging to the planar interface 

(see Fig. ||). Equation ( [T^ ) demonstrates that cr^gi°'^\d) 
never changes sign [as can also be deduced from simple 
geometrical considerations - Eq. (|l7|)] in contrast with 
the spherical interface. The total charge of polarization 

*3po/°"' obtained by direct integration of CTpo/""'' (d) and 
its expression is given by 



(sph) 
'pol 
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are the polarizations in region 1 and 2, respectively. Us- 
ing Eqs. (|), (i), (|), dll) and (|l|), the final expression 



of CTpo^'''' reads 
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The net charge of polarization Qp^i — 

J^^2na^a^^^i'^\cos9)d{cos9) is zero,0 meaning that 
there is no monopole contribution as it should be. 



Qiplan) 



£1 ' 



(20) 



where 



£1 - £2 
I : 

£1 + £2 



(21) 



is the unique image charge located at the mirror position 
of q (see Fig. |2]). This non-zero monopolar contribution 
for the planar interface involves a stronger and longer 
ranged self-image interaction. 



D. Application to colloidal systems 

So far we treat in a rather general manner the physics 
of a point charge near a spherical dielectric interface. We 
now would like to apply the above theory to colloidal 
systems. In the remaining of this paper we suppose that 
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FIG. 2: Model for a microion (g) near a planar interface 
(E) separating the two infinite media characterized by e\ and 
£2- The imaginary spherical dielectric of radius a is shown 
for geometrical comparison with the setup of Fig. |l| This 
is a two-dimensional representation of the three-dimensional 
system. 



region 1 corresponds to water, so that we take ei — 80 
corresponding to the water dielectric constant at room 
temperature. To characterize the low permittivity of the 
colloid we consider here £2 = 2 so that Ae = 78. The 
little ion carries a charge q = Ze where e stands for 
the elementary charge and Z for its valence, and has a 
diameter a. An important quantity is 



a 

2 



(22) 



being the center-center distance of closest approach be- 
tween the colloid and the microion q. 
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FIG. 3: Polar profile of the surface density of polarization 
charge o-^J'i'^\6) in units of a^J^ — for different radial 

distances b of the excess charge q with £1 = 80, £2 = 2 and 
a = 7.5(T. 



The numerical computation of Eq. ( p^ was performed 
using a cutoff ImcuL-.— 300 in the Legendre space leading 
to high accuracy.E3 The plot of o-p^^i'^\0) for a = 7.5<t 
and b/(T = 8, 9 and 10 can be found in Fig. ^. One can 
clearly observe that cr^^^''^ (0) is strongly inhomogeneous. 
For small 9, o-^^gi^\d) is highly positive (i.e., it carries the 



'pol 

same charge sign as q) and decreases abruptly. The an- 
gle 9* [given by Eq. (^] where tXp*-^'''' {9) changes sign is 
16.9°, 29.5° and 37.4° for b/a = 8, 9 and 10, respectively. 
In parallel, by increasing b one drastically decreases the 
magnitude as well as the inhomogeneity of cTpg^'''' {9) . Re- 
call that for 6/a ^ 1 we have o'p^fi'^\9) cos 9. 

It is insightful to compare <y^poi^\o = 0) with 
CTpg'""^ {d = 0) [computed from Eq. (p^] since both quan- 
tities give the maximum of Cp^^''^ {9) and cTpg;""^ {d) , re- 
spectively. The corresponding numerical values are gath- 
ered in Table |. The values found at finite curvature are 

isph) in _ 



very similar to those of zero one. The fact that a^^i 

0) is systematically smaller than cr^^j""^ {d = 0) is consis- 
tent with the idea that in spherical geometry we have 



1. Induced surface charge 

It is helpful to have a precise representation of the polar 
profile of CTpo^''^ {9) in order to get a clear understanding 
of the source of the image forces. Although at first glance 
such a study should belong to standard electrostatics we 
are not aware of any data in the literature that treats 
this crucial aspect. 



TABLE I; Numerical values of ol'^i\e = 0) and a^ir'id = 
0) in units of as a function of 6. The corresponding 



profiles of a<^„f) (e) can be found in Fig. H. 
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0.278 


0.304 
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the presence of opposite image charges. Nevertheless, for 
sufficiently large a one should recover the planar case. 



2. Self-image interaction 

We now compute the potential of interaction between 
the microion q and the dielectric particle or, in terms of 
image forces, the potential of self-image interaction. This 
is the work done in bringing the microion from infinity 
to its position b, and it is equal to the half -piodnct of q 
and the second term of $i(r = b) given by Eq. (0). In 
that case we have r = b (sec Fig. so that 9 = and 
therefore Pi[cos{9 = 0)] = 1. In order to normalize the 
energy with kgT we introduce the Bjerrum length Ig = 
eV(47reo£ifcBT) which is 7.14 A for water at T = 298 K. 
By choosing a = 3.57 A we have Ib = 2cr. The potential 
of self-image interaction V^^]'^'^ (b) is then given by 



v!:f\b) = hsTis^ 
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Equation ( p3| ) shows that the typical interactiop, range 
scales like 1/6"* and therefore it is short-rangedS3 Note 
that it is fully equivalent to compute V^^f^'^ (b) from the 
surface polarization charges as follows 
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where is the radial vector of magnitude |ra| = a and 

o-po^^ (cos 61) is given by Eq. (|l|). 

It is insightful to compare the potential of self-image 
interaction obtained in spherical geometry with that, 
vjfi^^\b — a), obtained in planar geometry. The setup 



for a planar interface is sketched in Fig. H. In this situa- 
tion the J 
given by 



tion the analytical expression of vj^i^"\b — a) is simply 
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FIG. 4: Potential of self-image interaction for a microion 
{q = Ze) in spherical and planar geometries with ei — 80 and 
£2 = 2. 
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self 



^^ij: (ro = 8a) = OMZ^keT and 
= 40cr) = OMZ'^ksT for the spherical inter- 

',^!iP\a/2) = 0.95Z^kBT for the planar one. 
These features can be physically explained in terms of po- 
larization charges. In the contact region (i. e., for small 
6 - see Fig. S) we know that the surface polarization 
charge is quasi-identical on both spherical and planar 
interfaces. However, for finite curvature we also know 
d{cos9), that CTp*^'*'' changes sign above 9* and in the present case 

(24) ^pli^^ sets oppositely charged to q. This latter effect is 
the main cause that leads to a weaker self-image interac- 
tion for spherical interfaces. Nevertheless, by increasing 
a (i.e., reducing the curvature) one approaches the pla- 
nar case as expected (see Fig. 0). Physically, this means 
that the contribution of the negative polarization charges 
(lying at 9 > 9*) to the self-image interaction [Eq. ([24|)] 
becomes negligible for sufficiently large colloidal radius. 



3. Effect of curvature on the contact potential 



1, rp, ^2^1^ £2 1 

— KrI IrZ ; -. 

2 " " £1+ £2 2(6 -a) 



(25) 



Profiles of V^'^i^\b) (for two colloidal radii) and 

^seif^^ (r) are reported in Fig. |4| Since in both (planar 
and spherical) cases the potential of interaction diverges 
at the interface, we only show results from r > tq as it 
is the case in experimental systems. The numerical com- 
putation of Eq. ( p3| ) was performed using thfi formalism 
of Eq. dl^) allowing an arbitrary precision.E3 

Figure h clearly shows that the self-image interac- 
tion is weaker (the higher the curvature) with a spheri- 
cal interface than with a planar one. In particular, at 



It is clear that for sufficiently low curvature one should 
recover the planar case as far as the self-image interaction 
is concerned. Thus, a natural question that arises is: for 
which typical colloidal size are curvature effects relevant? 

A suitable observable for this problem is provided by 
the contact potential vj^l'J^\b = a + a/2). This quan- 
tity is of special interest since it will correspond to the 
highest repulsive part of the global interaction between a 
macroion (i.e., charged macro-particle) and an oppositely 
charged counterion. In order to investigate the effect of 
finite curvature we are going to compare vj^fJ^\a + a/2) 

to the contact potential V^^ij^\b — a = cr/2) obtained 
with a planar interface. 
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The plot of the normahzed contact potential Vq* (a) de- 
fined as 



V*ia) 
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The plot of 
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can be found in Fig. |. For the sake of numerical stabil- 
ity we used the-£ormalism of Eq. (12) allowing an arbi- 
trary precision.c3 Figure ^ shows that for a/cr larger than 
about 100 the contact potential is close to that of the pla- 
nar interface (less than 5% difference). This length scale 
typically corresponds to "true" colloidal systems (~ 100 
nm). Therefore, in the dilute regime where the self-image 
interaction is dominant (i.e. lateral microion-microion 
correlations are negligible) , large-sized colloidal particles 
can be reasonably approximated by planar interfaces as 
far as the modeling of the self-image interaction is con- 
cerned. On the other hand, for a/ a smaller than about 
20 the contact potential varies rapidly and therefore it is 
strongly dependent on the curvature. This length scale 
typically corresponds to micellar systems (~ 10 nm). 

In this respect, Linse0 used an uncontrolled approx- 
imation where he replaced the (exact) infinite manifold 
of image charges [entering Eq. (|l2|)] of total charge —qim 
by a single image point-charge —<lirn-\ [given by Eq. (13)] 
located at the center of the sphere.Ell Doing so the setup 
of image charges consists of a (two point-charge) dipole 
Pim = Qim Uj and the corresponding contact potential 
Kil/^^(^o) reads 



Fj^^^)(a+f) 



(28) 



can also be found in Fig. |^. It shows that the two- 
image charge approximation used by Linse is only valid 
for very low curvature (i.e., close to the planar case) and 
may strongly overestimate the self-image repulsion as ex- 
pected by i1js-|inherent construction.Ej In his MC simula- 
tions, Linsell^ investigated micelles of radius 12 — 18 A 
(i.e, a/a ^ 3.5 — 5) leading to errors as large as 40% 
(see insert of Fig. This proves that this ansatz is 
unsuitable to determine the self-image interaction in this 
regime, which is the source of the image forces. 



4- Charged colloid 

As a last theoretical result, we consider the interaction 
between (a single counterion) q and a negatively charged 
dielectric sphere. The procedure is completely similar to 
the neutral colloid case, and we now apply the principle 
of superposition to take into account the additional po- 
tential due to a central charge Qm — ^ZmC. The (global) 
macroion-counterion potential of interaction (b) reads 




20 40 60 80 100 120 140 160 180 200 
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FIG. 5: Reduced contact potential Vo{a) as a function of 
the colloidal radius a with ei — 80, £2 = 2. The limit value 
of unity corresponds to the planar interface. The solid line is 
the exact contact potential Vo*(a) and the dashed one is the 
contact potential Vo'(t[) obtaiiped with the two-image charge 
approximation used by Linse t£l The insert shows the ratio 
Vo'{a)/Vo*{a). 
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where v}^l'!^\b) is given by Eq. (|2^), and hence 
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Profiles of Vm{b) for Zm — 60, ro — 8cr, £2 = 2 and 
Z = 1, 2 and 3 are reported in Fig. ^. An important re- 
sult is the occurrence of a minimum in Vm {b) whose depth 
and position r* increase with increasing Z. This is due 
to the purely repulsive self-image interaction which scales 
like Z^, whereas the direct attractive Coulomb macroion- 
microion interaction scales like Z (at fixed Zm)- Never- 
theless the occurrence of a minimum is strongly dictated 
by the ratio Z^/Z [see Eq. (|3C|)]. For high value of 
Zjn/Z^ \V^(b)\ is maximal for b — r^ (only attraction 
occurs) and for small Z^/Z one recovers the neutral col- 
loid case where only repulsion occurs. Of course the same 
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FIG. 6: Global macroion-counterion potential of interaction 
(solid lines) with Zm = 60, ro = 80", ei = 80 and £2 = 2. The 
values of the corresponding minima r* can be found in Table 
1^. The dashed lines correspond to the usual electrostatic 
potential of interaction without image forces (i.e., Ae = 0). 



qualitatively happens for charged plates .a The values of 
r* minimizing V^{\)) (with h > tq) are given in Table Q. 
The quantity r* will be useful to discuss our simulation 
results that concern many counterions and where we also 
have the same macroion bare charge {Zm = 60). 

Keep in mind that all our results above concern a single 
microion. When many counterions come into play, other 
important effects might appear in principle. In particu- 
lar, when the number of counterions near the macroion 
surface is very large the irnage forces are practically can- 
celed by symmetry reason.cj Clearly, by approaching the 
(perfect) spherical symmetry one asymptotically cancels 
the polarization charges everywhere on the macroion sur- 
face. This point shows that the discrete nature of the 
counterions is crucial_£or the existence of image charges 
in spherical geometry.E3 In planar geometry the situation 
is radically different, where one gets an amplified image 
force upon increasing the number of "surface" counteri- 
ons. 



III. MONTE CARLO SIMULATION 

Standard canonical MC aigmlations following the 
Metropolis scheme were usedrSa The system we^jGop. 
sider is similar to those studied in previous works. uO'L^I 
It is made up of two types of charged hard spheres: (i) 
a macroion of radius a with a bare charge — —Z^e 
(with Zm > 0) and (ii) small microions (counterions and 
colons) of diameter a with charge q — ±Ze to ensure the 
electroneutrality of the system. All these ions are con- 
fined in an impermeable cell of radius R and the macroion 
is held fixed at the center of the cell. 

The dielectric media are modeled as in Sec II. The 
outer region of the simulation cell is assumed to have 
the same dielectric constant ei as the solvent in order to 
avoid the appearance of artificial image forces. 

The work done in bringing the (real) ions together 
from infinite separation gives the interaction energy of 
the system. The corresponding Hamiltonian, Utot, can 
be expressed as 
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The first two terms in Eq. ( ^l| ) correspond to the tra- 
ditional electrostatic interactions between real charges. 
More explicitly. 
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represents the macroion-microion interaction, where (-I-) 
applies to colons and (-) to counterions, and 



±lBkBT — , for Tij > (T, 



00, 
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(33) 



the pair interaction between microions j and i where (-I-) 
applies to microions of the same type and (-) otherwise. 

The two last terms in Eq. (^ij) account for the in- 
teraction between images and microions. The repulsive 
self-image interaction is given by 



TABLE II: Theoretical values of r* minimizing the macroion- 
counterion potential of interaction (with Zm ~ 60, ei = 80, 
£2 = 2 and ro = Scr). The corresponding profiles can be found 
in Fig. |. 



£2 


A£ 


Z 


(r* - ro)/a 


2 


78 


1 





2 


78 


2 


0.17 


2 


78 


3 


0.32 



rrisel f) / \ 



21+1 



-kBTiB- y: 

' 1=1 

X — ^ — — , for > a + — , 
ei{l + l)+e2l 2' 



00, 



for ri < a + 



(34) 
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where 



is the cutoff in the Legendre space, and 



there is no dielectric discontinuity (Ae = 0). Measure- 
ments were performed over 10^ MC steps per particle. 



-,21+1 



{£l - 62)1 



"hi £i(/ + l) + £2Z 

x-j-jPi{cos9), forri>a + ^, 

2, 

a 
2' 



oo, 



for Ti < a 



(35) 



represents the interaction between microion i and the 
image (surface charge induced by) of microion j, where 
(+) applies to charges of the same sign [and (-) otherwise] 
and 9 is the angle between r j and Vj . It is this term that 
generates lateral image-counterion correlations. Due to 
the symmetry of C/j^™^ upon exchanging ij with ji there 



is an implicit factor 1/2 in Eq. (pS]). 

Convergence of the Legendre sums with a relative 
error of 10~p£, is obtained with the employed value of 
Imax = 100.t3 For the sake of computational efficiency 
and without loss of accuracy, we computed the image-ion 
interactions on a (very) fine (r, cos 9) grid where the co- 
ordinates of the microions were extrapolated. The radial 
distances are discretized over logarithmically equidis- 
tant nodes so that close to the macroion surface the radial 
resolution is O.Olcr and near the simulation wall O.la. The 
polar discretization consists of 2000 equidistant cos 61- 
nodes leading to even smaller lateral resolutions. The 
corresponding values of [/.^"'^'•^^ (r^) and U^^™^ {ri , rj , cos 9) 
were then initially stored into tables. Note that it in prin- 
ciple one could also have used the formalism of Eq. ( [l^ ) 
to compute the image-ion interactions. However, at iden- 
tical numerical accuracy, this method involving a numer- 
ical integration is too time and resource consuming. 

Typical simulation parameters are gathered in Table 
[II. The case £2 — 80 corresponds to the situation where 



TABLE III: Model simulation parameters with some fixed 
values. Apart from the charge sign, counterions and coions 
have the same parameters. 

Parameters 



T = 298is: 
£1 = 80 

£2=2 

A£ = £1 - £2 = 78 

z 

a = 3.57 A 
lB=2a = 7.14 A 
a = 7.5(7 

ro — a + ^ — 

R 



room temperature 
water solvent dielectric constant 
colloidal dielectric constant 
strength of dielectric discontinuity 
macroion valence 
counterion valence 
counterion diameter 
Bjerrum length 
macroion radius 
macroion-counterion distance 
of closest approach 
radius of the outer simulation cell 



IV. SIMULATION RESULTS 

Here we present our MC simulation results in salt-free 
environment as well as in the presence of multivalent salt- 
ions. We essentially study in detail the radial microion 
distributions ni{r) around the macroion, which are nor- 
malized as follows 



47rr^n+(r)(ir — N-\- 



(36) 



where r is the distance separation from the macroion cen- 
ter, -l-(-) stands for counterion (colon) species and 
(iV_) is the total number of counterions (coions) con- 
tained in the simulation cell. 

Another quantity of special interest is the integrated 
(or cumulative) fluid net charge Q{r) defined as 



Qir) 



Attu'^Z [n+{u) - n_(w)] du, (37) 



where we chose e — 1. Q{r) corresponds to the total fluid 
charge (omitting the macroion bare charge Z^) within a 
distance r from the macroion center, and at the cell wall 
Q{r = R) ^ Zjn- Up to a factor proportional to l/r^, 
[Q{r) — Zm] gives (by simple application of the Gauss 
theorem) the mean electric field at r. Therefore Q{r) can 
measure the strength of the macroion charge screening 
by salt-ions. In salt-free environment systems we have 
n_(r) = and N+ = Z^/Z. 

The simulation run parameters can be found in Ta- 
ble For all these simulation systems, the ion den- 
sities riiW were computed with the same radial resolu- 
tion Ar.cJ The discretization of the radial distance r in 
ni{r) is realized over logarithmically equidistant points 
so that close to the macroion surface (r — ro < a) we 
have Ar < Q.OAa. It is important to obtain such an 



TABLE IV: System parameters. 



System 


A 


B 


C 


D 


E 


F 


G 


H 


/ 


J 




60 


60 


60 


60 


60 


60 


60 


60 


180 


180 


Z 


1 


1 


2 


2 


3 


3 


2 


2 


2 


2 


N+ 


60 


60 


30 


30 


20 


20 


430 


430 


445 


445 


N- 














400 


400 


400 


400 


£2 


2 


80 


2 


80 


2 


80 


2 


80 


2 


80 


A£ 


78 





78 





78 





78 





78 





R/a 


40 


40 


40 


40 


40 


40 


20 


20 


20 


20 
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accuracy (and the required statistics) if one wants to de- 
scribe quantitatively the effects of image forces which are 
short-ranged at strong curvature. 



Salt-free environment 



0.04 



0.03 




Sah-free systems A — F (see Table fV) were investi- 
gated for a moderately charged macroion = 60 cor- 
responding to a surface charge density ctq — 0.11 Cm~^. 



s 



0.02 



0.01 



Monovalent counterions 



The profiles of n^{r) and Q{r) are depicted in Fig. 
0(a) and (b), respectively for the monovalent counterion 
systems A and B. 

Figure ^a) shows that the counterion density at con- 
tact (r = To) is somewhat smaller with Ae = 78 as a 
direct consequence of the self-image repulsion. However 
there is no maximum appearing in n^{r) with Ae = 78, 
in agreement with the study of the single-countcrion sys- 
tem (see Fig. I and Table H). For r — ro >~ 0.6a 
(corresponding roughly to three half ionic sizes from the 
interface), the effects of image forces are negligible and 
all n+(r) curves are nearly identical. 

To gain further insight into the effects of lateral image- 
counterion correlations, we have considered the same sys- 
tem A {Ae = 78) but omitted the correlational term 

C7^™^ [Eq. (H)] in the total Hamiltonian f/tot [Eq. (|l|)]. 
Physically, this means that, on the level of the image 
force, each counterion sees uniquely its self-image interac- 
tion. Thereby, Fig. ^(a) shows that (i) the corresponding 
counterion density rS^^^^^{r) is nearly identical to n+(r), 

and (ii) in the vicinity of the interface ■n}-_^'^^^'^ (r) is slightly 
smaller than n+ (r) . These findings (i) and (ii) lead to the 
two important conclusions: 

• For monovalent counterions and moderately 
charged macroions, the effective image force 
is basicallyj— jdentical to that of the self-image 
interaction.LJ 




2 3 

(r-r„)/a 

Z = 60 ; Z = 1 




FIG. 7: Monovalent counterion distributions (systems A and 
B): (a) density n+(r). The dashed line in grey corresponds 
to the counterion density n^^'^^^\r) obtained in the same sys- 
tem A (Ae = 78) but where the (lateral) image-counterion 



(zm) 



correlational term 
total Hamiltonian Utot [Eq. 



[Eq 



] has been omitted in the 
(b) fluid charge. 



• The crucial effect of lateral image-counterion cor- 
relations is to screen the self-image repulsion. 

This latter feature is generally true for any finite curva- 
ture at identical fixed macroion charge density. Finding 
(i) is also consistent with the fact that, close to the inter- 
face (say r — ro < 0.2(t), the average number of (surface) 
counterions N is (very) small (A^ < 5) as can be deduced 
from the fraction of counterions Q{r)/Zm [Fig. ^b)]. 

Figure 0(b) shows that the fluid charge Q{r) decreases 
when image forces are present, meaning that they lower 
the macroion charge screening by counterions. At the 
distance r — rp — a (corresponding to a 2cr-layer thick- 
ness), the macroion is 29% electrically compensated [i.e., 
Q{r-ro = a)/Zm = 0.29] with Ae = against 26% with 



Ae = 78. At the distance ? — tq = 4a, the relative differ- 
ence AQ/Q between the Q{r) obtained with Ae = and 
Ae = 78 drops to 2% (against 10% at t — ro = a) where 
the bare macroion charge is nearly half-compensated. 



2. Multivalent counterions 

a. Divalent counterions The profiles of ?i+(r) and 
Q(r) are depicted in Fig. ^(a) and (b), respectively for 
the divalent counterion systems C and D. 

Figure ||(a) shows that the counterion density at con- 
tact becomes strongly reduced with Ae = 78 due to the 
Z^-dependence of the self-image repulsion [compare the 
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Z = 60 ; Z = 2 



0.06 



0.05 



0.04 
^ 0.03 
0.02 



0.01 



Ae = 

Ae = 78 

Ae = 78 



12 3 4 

(r-r„)/a 



Z = 60 ; Z = 2 




FIG. 8: Divalent counterion distributions (systems C and 
D): (a) density n+(r). The dashed hne in grey corresponds 
to the counterion density n^,f'^''^'(r) obtained in the same sys- 
tem C (Ae = 78) but where the (lateral) image-counterion 



correlational term 
total Hamiltonian Utot [Eq. (pll)] 



[Eq 



] has been omitted in the 
(b) fluid charge. 



case Z = 1 vcL Fig. 0(a)]. This sufficiently strong (short- 
ranged) repulsion leads to a maximum in n-\-(r) close to 
the macroion surface. The corresponding radial position 
r* maximizing 7i+(r) is r* = -\- 0.22a, in excellent 
agreement (within Ar) with the one- counterion theoret- 
ical value To + 0.17(7 (see Table Q). This shows that 
for divalent counterions many-body effects do nearly not 
affect r* . This non-trivial finding is the result of the 
competition between two driving forces that control r* 
in manj/-counterion systems: 

• Fim- the screening of the self-image repulsion by 
the (extra) negative polarization charges tends to 
decrease the r* obtained in the one-counterion sys- 
tem. 



• Fmc- the screening of the macroion-counterion at- 
traction by the (extra) surface counterions tends 
to increase the r* obtained in the one-counterion 
system. 

It is precisely a balance of these two driving forces 
that leads to a nearly unchanged r* (compared to 
the one-counterion system) in many-counterion systems. 
Whereas for monovalent counterions both driving forces 
and Fmc are weak, those become relevant for multi- 
valent counterions. 

We stress the fact that this is specific to the spherical 
geometry, and that for a planar interface (at identical 
surface charge density) one should get a higher r* (com- 
pared to that of the one-counterion system), since there 
we have no screening driving force Fim ■ We are not aware 
of any previouSj-studies for the planar interface that ad- 
dress this issue .Lj 

To gain even further insight into the effect of Z on the 
lateral image-counterion correlations, we have ignored 
the term C/j^^""^ in Utot in the same system D (Ae — 78) 
as done previously with system A. Figure ^(a) shows a 
qualitatively different n";^'^^^^ (r) where r* = r^ -\- 0.26ct is 
now somewhat larger, proving that with divalent counte- 
rions the screening of the self-image repulsion by lateral 
image-counterion correlations is appreciable. This is in 
contrast to what was observed with Z = \. 

At the distance r — ro — a, Fig. ||(b) shows that the 
macroion is 62% electrically compensated for Ae = 
against 53% for Ae = 78 [compare the case Z — 1 in Fig. 
0(b)]. 

h. Trivalent counterions The profiles of n^{r) and 
Q{r) are depicted in Fig. ^(a) and (b), respectively for 
trivalent counterion systems E and F. 

Figure 0(a) shows that the counterion density at con- 
tact is drastically reduced with Ae = 78, as expected for 
high Z (compare the previous cases). At Ae = 78, we 
have r* = rg -I- 0.36(7, in quantitative agreement with the 
one-counterion theoretical value ro -I- 0.32(7 (see Table 
0). This shows again that even for trivalent counterions 
many-body effects do (practically) not affect r* (com- 
pared to that obtained in the single-counterion system) 
due to a balance of the driving forces Fim and Fmc- 

By neglecting the lateral image-counterion correlations 
in the same system E (Ae = 78), Fig. 0(a) indicates 
that the position r* of the maximum in n^^'^^^^ (r) gets 
considerably larger (r* = tq + O.SOcr). This relatively 
strong shift confirms the Z-cnhancing of the screening 
of the self-image repulsion by lateral image-counterion 
correlations. 

At the distance r — ro — a, the macroion is 84% 
electrically compensated for Ae = against only 67% 
for Ae = 78 [see Fig. 0(b) and compare previous sys- 
tems]. Snapshots of typical equilibrium configurations 
for Ae = and Ae = 78 can be visualized in Fig. |l^ (a) 
and (b), respectively. 
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FIG. 10: Snapshots of typical equilibrium configurations 
for trivalent counterions (systems E and F). (a) Ae = 
(b) Ae = 78. One can clearly observe the larger mean radial 
counterion distance for Ae = 78 stemming from the self-image 
repulsion. 



1. Moderately charged macroion 



FIG. 9: Trivalent counterion distribution (systems E and 
F): (a) density n+(r). The dashed line in grey corresponds 
to the counterion density ri-^'^^^\r) obtained in the same sys- 
tem E (Ae = 78) but where the (lateral) image-counterion 
correlational term C/^'j™' [Eg. ( ] 35[ )] has been omitted in the 
total Hamiltonian Utot [Eq. (|31|)]. (b) fluid charge. 



B. Salty solutions 

We focus on the case of divalent salt-ions. This choice 
is motivated by two reasons: (i) effects of image charges 
are clearly observable for multivalent counterions and (ii) 
such systems must be experimentally reachable. To study 
the effect of added salt we have considered two macroion 
charges Zm = 60 (as previously) and Z^n — 180 corre- 
sponding to a charge density ctq — 0.32 Cm^^. The salt 
concentration defined as ^^"3 is 0.44 M for all salty sys- 

tems G — J (see Table [IV| ). The simulation cell radius 
R = 20(J of these systems is still very large compared to 
any screening lengths so that finite size effects are negli- 
gible. 



Profiles of n±{r) and Q{r) are depicted in Fig. [Tl](a) 
and (b) , respectively for the salty systems G and H with 

Z,n = 60. 

The colon density n_(r) with Ae = 78 is basically 
shifted to the right of about O.lScr (compared to that 
with Ae — 0) due to the repulsive coion' self-image in- 
teraction. Near the colloidal surface, the counterion den- 
sities n^{r) are considerably higher than those obtained 
with no added salt (systems C and D) as it should be 
[compare Fig. ^(a)]. 

A rather surprising result here is that, despite of the 
presence of a considerable amount of added salt, we still 
have r* = + 0.22a remaining unchanged. This is a 
non-trivial finding since one should have an (extra) at- 
tractive contribution to the macroion-counterion poten- 
tial of mean force stemming from the (localized) negative 
polarization charges induced by the colons, which in turn 
could lead to a shorter r* . However there are two con- 
comitant sources that lead to a marginal screening of the 
counterion' self-image repulsion by the negative coion- 
induced polarization charges: (i) there is a strong coion 
depletion close to the interface [see Fig. ^l](a)] due to 
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Z = 60 ; Z = 2 



Z = 180;Z=2 



0.08 



0.06 



^ 0.04 
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(r-r„)/a 
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Ae = 



Ae = 78 

Ae = 78 
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FIG. 11: Divalent salt-ion distribution (systems G and H) 
with Zm = 60: (a) The solid and dashed lines correspond 
to counterion and colon densities, respectively, (b) net fluid 
charge. 



FIG. 12: Divalent salt-ion distribution (systems I and J) 
with Zm = 180: (a) The solid and dashed lines correspond 
to counterion and coion densities, respectively, (b) net fluid 
charge. 



the large direct Coulomb macroion-coion repulsion and 
(ii) |o'po^''^| decreases abruptly with the radial distance of 
the microion as discussed in Sec. [ID 1 (see also Fig. ||). 
Of course the role of the excluded volume is crucial here. 

As expected the macroion charge screening is weaker 
when image forces come into play as can be deduced from 
the profile of Q{r) plotted in Fig. |ll|(b). 



2. Highly charged macroion 

Profiles of n±{r) and Q{r) are depicted in Fig. [T^(a) 
and (b), respectively for the salty systems / and J with 

Zm = 180. 

Figure |lj(a) shows that the effects of image forces are 
considerably reduced. The relatively small difference be- 



tween the n+ (r) obtained with Ae — and that obtained 
with Ae = 78 decreases drastically in the vicinity of the 
interface, and already for r — vq 0.2a the two profiles 
of n+{r) are nearly identical. Besides, near the inter- 
face no effective macroion-counterion repulsion occurs at 
Ae = 78. This absence of a maximum in n^{r) is due to 
two main concomitant effects: 

• For such a highly charged macroion, there is a very 
large number of counterions close to the interface 
[compare Fig. |l|(b) and Fig. 0(b)]. In this limit, 
one can use Wigner crystal concepts and say that, 
on the level of the force stemming from the bare 
charges (i.e, ignoring the image forces), each surface 
counterion essentially interacts with the oppositely 
charged background of its Wigner-Seitz (WS) cell. 
At sufficiently high macroion charge density (i.e.. 
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small WS hole radius), this attractive interaction 
becomes very important and it always overcomes 
the self-image repulsion. 

• The second (concomitant) mechanism is specific to 
the closed spherical topology: at high number of 
surface counterions, the image forces are reduced 
because of the enhanced degree of spheri cal sym- 
metry as already mentioned in Sec. [ID 4 



The colon densities n_(r) are basically identical for 
both dielectric discontinuities As, in contrast to what 
happened with Zm = 60 (systems J and K). This non- 
trivial finding can be explained as the enhanced screening 
of the colon' self-image repulsion by the positive polar- 
ization charges induced by the other colons present in 
the electrical double layer (EDL). Indeed, because of the 
macroion charge reversal that occurs at Zm = 180 [i.e., 
Q{r)/Z„i > 1 - see Fig. |l^(b)], there is also a larger 
number of colons (at fixed salt concentration) in the EDL 
[compare Fig. |l|(a) and Fig. |ll|(a)]. Therefore, since the 
magnitude and the inhomogeneity of — cr^^^'''' [0) induced 
by a colon strongly decreases with its radial distance [see 
Eq. ( |l^ ) and Fig. |j, the screening of the colon' self- 
image repulsion gets highly sensitive to an increase in 
number of colons in the EDL. 

Concerning the net fluid charge (5(r), we see that both 
profiles obtained with Ae = 78 and Ae = are nearly 
identical, as expected from those of n± (r) . The net fluid 
charge Q{r) reaches its maximum at Tq—r^ = 0.90ct and 
0.94a for Ae2 = and 78, respectively. In both cases 
we have a macroion charge reversal of 9% [more explic- 
itly Q{rQ)/Zm — 1.09]. This proves the important result 
that, for typical systems (with high mijrrpion charge den- 
sity) leading to overchargingJ^OcJ'EllEj image forces do 
not affect the strength of the macroion charge reversal. 



V. CONCLUDING REMARKS 

We have presented fundamental results about the 
effects of image forces on the counterion distribution 
around a spherical macroion. 

Exact analytical results have been provided for the 
case of a single microion interacting with a dielectric 
sphere. Within this framework, the self-image interac- 
tion and the surface charge of polarization have been 
studied and also compared to those obtained with a pla- 
nar interface. Besides we also estimated the position r* 
where the macroion-counterion potential of interaction is 
minimized. We demonstrated that the effects of image 
forces due to a spherical interface are qualitatively dif- 
ferent from those occurring with a planar interface, espe- 
cially when the colloidal curvature is large. We showed 
that the self-screening of the polarization charges (i.e., 
the screening of the positive surface charges of polar- 
ization by the negative ones) is decisive to explain the 
weaker and the shorter range of the self-image interac- 



tion in spherical geometry. This self-screening increases 
with the colloidal curvature. 

Many-counterion systems have been investigated by 
means of extensive MC simulations where image forces 
were properly taken into account. 

In salt-free environment and for moderately charged 
macroions, a maximum in the counterion density (near 
the spherical interface) appears for sufficiently large di- 
electric discontinuity Ae. An important result is that the 
corresponding position r* is basically identical, regardless 
of the counterion valence Z , to that obtained within the 
one-counterion system. This feature is specific to the 
spherical geometry and can not take place with planar 
interfaces where there is no self-screening of the polariza- 
tion charges. For monovalent counterions we showed that 
the (effective) image force is basically equal to that of the 
self-image interaction, and the lateral image -counterion 
correlations are (very) weak. However for multivalent 
counterions the lateral image-counterion correlations af- 
fect significantly the counterion density, and as major 
effect they screen the self-image repulsion. Nevertheless, 
the combined effects of (i) the macroion charge screening 
by counterions and (ii) the screening of the self-image re- 
pulsion lead to a nearly unchanged r* (compared to that 
obtained in the single-counterion system) for multivalent 
many-counterion systems. We also showed that the coun- 
terion density at contact decreases drastically with Z, 
and that r* also increases with Z as expected. These 
latter results have important implications for the stabi- 
lization of charged colloidal suspensions where a compo- 
nent of the pair-force is proportional to the ion density 
at contact. 

By adding salt, it was found for moderately charged 
macroions that the strength of the image forces induced 
by the coions is very small compared to that resulting 
from the counterions. This is due to the coupled effects 
of (i) the colon depletion in the vicinity of the colloidal in- 
terface due to the strong direct Coulomb macroion-coion 
repulsion and (ii) the (highly) short range of the image 
forces in spherical geometry. Consequently the position 
r* remains identical to that obtained in salt free envi- 
ronment and a fortiori to that obtained within the one- 
counterion system. For highly charged macroions the ef- 
fects of image charges are significantly reduced since (1) 
the attractive counterion-hole interaction dominates the 
repulsive counterion' self-image interaction and (ii) the 
screening of the counterion' self-image repulsion gets en- 
hanced by symmetry reason. In this situation no maxi- 
mum appears in the counterion density and it was found 
that overcharging is nearly unaffected by image forces. 

Although our MC analysis was carried at given 
macroion size, all the above reasonings that concern 
many counterions remain unchanged (for symmetry rea- 
son) for any finite curvature by a rescaling at fixed 
macroion charge density. 

Finally, this contribution should constitute a solid ba- 
sis to understand and predict the effects of image charges 
in other similar systems (e.g., polyelectrolyte adsorption 
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onto spherical charged colloids). 
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thogonality of the Legendre polynomials. Indeed, we have 
/_i Pi(x)Po{x)dx = and since Po{x) = 1 it results from 

this that j^^Pi{x)dx = 0. 
^® It turns out that the convergence of the Legendre sum 
in Eq. ( [l^ ) is extremely slow. With Imax ~ 300 we have a 
relative error of 10~^ in the most severe case where b = 80". 
This holds when the microion q is sufficiently far away 
from the interface so that only the dipolar term matters 
(i.e, / = 1). However close to the colloidal surface one has 
to take all higher multipolar contributions to get the exact 
potential of interaction. 

We have carefully checked that both methods [i.e., Eq. ( px| ) 
and Eq. {^] are numerically identical. However, at large a 
(typically a > 30a) the Legendre sum in Eq. (ESI) becomes 
numerically unstie^le. 

Originally, Linsec3 made some manipulations in the Leg- 
endre expansion and drops all I components after having 



alents to replace the ("exact) infinite manifold of image 
charges [entering Eq. (jl^)] by a single image point-charge 
—qim located at the center of the sphere. 
Note that the two-image charge approximation also holds 
when the charge q is very far from the interface (i.e., b/a ^ 
1). Unfortunately, all the relevant effects of image forces in 
spherical geometry occur within the opposite limit b/a 
1. 

In the case of a charged plate, to be consistent, one has 
to think either of a large finite plate in the presence of a 
single microion, or equivalently of an infinite plate in the 
presence of a very dilute dispersion of microions. 
C. Holm (private communication). 

This problem is very similar, in spirit, to that of the dis- 
cretization of the macroion surface charge whose effect is 
(asymptotienkll*') canceled for small separations between 
microions. EflEZi 
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Here, the convergence of the Legendre sum in Eqs. (^) 
and ( p5[ ) is much faster, in contrast to what happened with 
the computation of o"po?'''(^) [Eq. ([l6|)]. The most sensitive 
region is the one where r — » ro. There, for ro = 8a we have 
a relative error of 10~® with Imax = 100. 
For the net fiuid charge Q{r) one has an arbitrary resolu- 
tion Ar. 

This feature also shows how important it is to treat pre- 
cisely the self-image interaction as done in the present 
study with a Legendre sum. n 

Using MC simulations, Torrie et ala also reported in the 
planar case the appearance of a maximum in the counte- 
rion density for low surface charge density and fairly small 
salt content. However, they did not make any comparison 
between the position of this maximum and the position 
of a single counterion (per simulation box) minimizing the 
plate-counterion potential of interaction. 



